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Qh. ABSTRACT: We establish the phase diagram of three-dimensional quantum gravity 

^ ■ coupled to Ising matter. We find that in the negative curvature phase of the quantum 
gravity there is no disordered phase for ferromagnetic Ising matter because the coordination 
^ ! number of the sites diverges. In the positive curvature phase of the quantum gravity there 
is evidence for two spin phases with a first order transition between them. 

One of the results that propelled recent work on two dimensional quantum gravity 
was the discovery that dynamically triangulated random surfaces provide a convenient 
formulation of the theory and facilitate both analytic and numerical calculations. An 
analogous formulation of three dimensional quantum gravity has also been proposed [1] 
with some success [2]. This success is non-trivial because the proof for two dimensions 
that the number of triangulations is exponentially bounded as a function of the number of 
triangles does not carry over to three dimensions. In two dimensions, the classification of 
topologies in terms of genus is crucial in proving the bound. In three dimensions, there is 
no comparable topological classification scheme. Without the exponential bound, the naive 
result is to have factorial growth. This would mean that no chemical potential term could 



control the theory. However, numerical experiments suggest that there is an exponential 
bound, at least when the topology is restricted to [3] . 

Another difference in three dimensions is the need for a second chemical potential. In 
two dimensions, the Euler characteristic, Xi is 

No-Ni + N2 = x (1) 

(where Nj^ is the number of i-simplices). Furthermore, each triangle has three links and 
each link is shared by two triangles implying 

3N2 = 2Ni (2) 

The most general action of the form 

S = Y,9iNi (3) 

therefore has only one independent term. This can be taken to be a cosmological constant 
determining the mean number of triangles. In two dimensions it is further possible to 
restrict ourselves to a "microcanonical" ensemble which eliminates even this term. The 
corresponding equations in three dimensions are (for x — which will be assumed from 
here on) 

iVo - iVi + iV2 - A^3 = (4) 

and 

4iV3 = 2N2 (5) 

leaving two independent terms in an action of the form (3). This action can be chosen to 
be 

S = aNo - /3N'i (6) 

so that (3 controls the number of tetrahedra (the volume) and a controls the number of 
nodes. In the continuum, /3 corresponds to the cosmological constant while a is related 
to the gravitational constant and hence the mean curvature. Unfortunately, there are no 
proofs of ergodicity for a microcanonical ensemble in three dimensions, so we are forced to 
utilize a canonical ensemble for our simulations. Of course, the existence of a continuum 
limit is a major issue. The question is whether or not there is a second order phase 
transition in the (a, /?) plane (or in some larger space if actions of a more general form 
than (3) are considered). The coupling /? must be adjusted so that the volume approaches 
infinity. This leaves a for our attention. It is reported that for a < a* (where a* ~ 
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3.75 — 4.0), the theory is in a negative curvature phase while for a> a* the mean curvature 
is positive [4] . Simulations suggest that the transition between them may be first order [5] 
implying that no continuum limit is available. It has been speculated that this negative 
result is related to the topological nature of pure three dimensional gravity in which there 
is no massless graviton. 

Matter can be coupled to the gravity by placing it at the nodes (for instance) and 
adding an appropriate term to the action. In two dimensions, there is some understanding 
of the effects of coupling quantum gravity to matter, at least when the central charge is 
less than one. Knihznik, Polyakov, and Zamalodchikov (KPZ) [6] have argued that when 
the central charge is less than one, the scaling dimensions of the matter fields are simply 
"dressed" by quantum gravity. On the lattice this means that for a critical point of some 
matter theory on a fixed lattice, there is a corresponding second order phase transition 
for the matter on a dynamical lattice, but with different critical exponents. KPZ give an 
explicit relation for the two sets of scaling dimensions. In particular, the Ising model has 
a continuous phase transition when coupled to quantum gravity with exponents different 
from the usual ones. In this paper, we couple the Ising model to three dimensional quantum 
gravity [7] . It is possible that within this enlarged parameter space we may be able to find 
new critical regions where a contimmm limit could be taken. There is no analog of the 
KPZ framework to guide us, so, clearly, a central issue is whether or not an Ising transition 
continues to exist and, if so, what its order is. 

We introduce the matter by simply decorating the nodes of each diagram with spins, 
adding an Ising term to the action with coupling A;, and tracing over spins ti = ±1. Thus 
our total partition function Z is given by 



There are two important relations between the partition function of this new theory 
and the partition function of pure gravity. The first is due to the fact that the trace over 
spins contributes a factor 2^° even if the Ising coupling is zero: 



In other words, in a grand canonical ensemble, the Ising spins renormalize the gravitational 
couplings even at zero coupling. For the new action, this implies that the curvature 
transition is in the range a* ~ 3.06 — 3.31. The second relation is due to the fact that 





Z{a, P,k = 0) = Zpure grav{a + In 2, (3) 



(8) 
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once the spins are completely ordered, the Ising term just contributes a factor A;A^i to the 
action. 



Zia, (5, k) ~ Zpure grav{a + k,P - k) (for k such that < tt' >r^ Ni) (9) 

One implication of this last result is that starting with any a in the negative curvature 
phase, large enough k will take the system into the positive curvature phase. In the 
discussion that follows we always assume that (3 has been adjusted to fix the mean volume 
associated with the tetrahedra. The continuum limit will finally be taken by letting this 
mean value tend to infinity. This leaves us with a two dimensional phase diagram. (This 
can be done because the empirical exponential bound for the pure gravity implies an 
exponential bound in the presence of matter). 

Before proceeding with a numerical analysis of the phase diagram it is useful to consider 
the Ising model in mean field theory on quenched lattices of the type found in the pure 
gravity theory. This will give a kind of zeroth order approximation to the phase diagram 
and will help direct the numerical work. 

Take A^o ^ind to be fixed and make the standard mean field approximation 

J2tt' = q<t>J2t (10) 

where tt' is the Ising term, < t > is a spontaneous magnetization, and q is the average 
coordination number. As usual we may add an external magnetic field to facilitate the 
computation of expectation values. The partition function now reads 

Z = (const) Ti exp{[h + kq < t >]^t) (11) 

where h is the magnetic field and k is the Ising coupling. This is now a non-interacting 
system so the trace can be done to get 

Z = (const) (2 cosh[/t + kq < t >])^° (12) 

The average spin can be calculated from this and is 

1 dlnZ , /, , s 

<t>= — -— = tanh(/i + kq<t>) (13) 

Nq ah 

For h = and small < t > (i.e. near the critical point) < t >= kq < t >, so the critical 
value of the Ising coupling is 

k = 1/q. (14) 
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Equations (10-14) are standard. What is new is the effect of the background lattice 
via the coordination number q. Each hnk contributes to the coordination number of two 
sites, so 

2Ni = Noq (15) 
Equations (14) and (15), together with A^i = Nq + (from (4) and (5)), give 

^ = \{l + {N,/No)) ^^^^ 

Now the problem boils down to the ratio of N^/Nq as A^3 goes to infinity. Ref. [4] finds 
that for large A^3, < Nq >= N^^^ with 5 < 1 in the negative curvature phase and 6 = 1 
in the positive curvature phase. When 5 < 1, N^/Nq goes to infinity (which just means 
the average coordination number goes to infinity) and the critical coupling is at zero. 
This is what happens at all points in the negative curvature phase. When 5 = 1, i.e. 
< A^o >= -^s/c, we get 

Ref. [4] gets (in pure gravity) for a = 4, c = 3.368. This gives /c = 0.114 (which turns out 
to be very close to our numerical result). 

This mean field argument suggests the phase diagram of fig. 1. The dashed line is 
just eqn. (9) with the knowledge that the negative curvature phase orders the spins for 
arbitrarily small couplings. This line then separates the negative curvature phase from the 
positive curvature phase. The thick line is equations (16) and (17) with input about the 
relationship between Nq and A^3 in the two phases. Since mean field theory is an expansion 
in 1/q, when g = oo it should be exact. In the negative curvature phase, therefore, the 
Ising spins order for arbitrarily small positive values of the Ising coupling. The only part 
of the phase diagram in question is the positive curvature phase part of the Ising line. We 
address this question numerically. 

The gravity sector of the theory is simulated using update moves that have been shown 
to be ergodic [8] and have become standard. Barycentric subdivision, also called the (1,4) 
move, takes a tetrahedra, inserts a node in the middle and connects this node to the 
other four vertices creating a total of four tetrahedra. The analog of the bond flip in two- 
dimensions, here the (2,3) move, removes the triangle shared by two tetrahedra, introduces 
a link between the two nodes that are not shared and places three new triangles around 
this new link (one for each of the shared nodes) in such a way that there are now three 
tetrahedra. 



There is a fine tuning problem associated with this theory. If the number of simphcial 
manifolds regarded as a function of the number of tetrahedra behaves as 

#(Ar3) ~ (Ar3)«e^^3 (18) 

with a negative (such as happens for two dimensions), then the grand canonical partition 
function is dominated by large volumes if the chemical potential is slightly too small and 
by small volumes if the chemical potential is slightly too large. In numerical simulations 
this manifests itself as the tendency for the system to evolve to one of the two extremes no 
matter how carefully the couplings are adjusted. Several approaches to this difficulty have 
been used. One group confined their simulations to a range of A^a [4]. We prefer to retain 
an action formulation of the problem and add a gaussian term to the action as proposed 
in [2] (except that we do not constrain A^o)- 

5S = -7(iV3 - Ma)^ (19) 

The desired volume M3 is put in by hand. The coupling j3 is then adjusted so that 
< >— M3. The new coupling 7 is made small so that it does not influence the physics 
of interest. This new term in the action effectively ensures that the probability distribution 

for A^3 is approximately gaussian about the mean value M3, with a characteristic width 
determined by 7. The length of the simulation is then tuned with 7 so that the system has 
time to wander back and forth across the distribution many times. Otherwise, the system 
may not reach equilibrium. In practice, 7 = 0.005 is a practical value. 

With the moves and action now defined, it is important that detailed balance must be 
built into the update probabilities [9]. For pure gravity, two lists are made, one with all 
sets of two or three tetrahedra for which the 23 or 32 move is allowed and another with 
all sets of one or four tetrahedra for which the 14 or 41 move is allowed. Each list consists 
of items in one of two possible states which we will call A and B. There is some number, 
n, of items in the list. Picking a state at random and updating it produces a modified list 
with n' items. Detailed balance is maintained if the updating probability is 

P{A ^B)= \ (20) 

Adding the Ising matter requires some modifications. Except for its effect on the action, 
no modification is necessary for the 23 or 32 moves. However, the update probabilities 
must be modified for the 14 and 41 moves because the list now has three possible states 

6 



(no node, node spin up, node spin down) rather than two. Call the spin states Bl and 
B2. A suitable modification of the update probabilities is 

P{A ^ Bl) 



P{B1 -^A) = ^ (21) 



P{B1 B2) 



ne^B2 



and three more expressions with 1 and 2 interchanged. It is most convenient to consider the 
total probability of introducing a new node, P{A Bl) + P{A — ^ B2), and to determine 
the spin value (using Metropolis) after deciding whether or not to put in the new node. 

Of course, there are other ways to maintain detailed balance. One can choose from lists 
of moves which are not necessarily allowed and then test for this during the accept/reject 
step. In the approach above, it is most time consuming to predict the number of elements 
that will be in the list of allowed (23) moves if a contemplated (23) move is made. It is 
more efficient to postpone part of this calculation to the accept/reject step (the part that 
checks whether or not a candidate simplex already exists or not) while keeping the simpler 
restrictions on the list (such as requiring links to have the correct coordination number). 

Calculations were performed for values of a in the range 3.5 to 4.0 and for values of 
the Ising coupling between and 0.2. As A; is varied in this range the nearest neighbor 
spin expectation value varies from to around 0.7 (i.e. from complete disorder to fairly 

ordered). The Ising component of the specific heat was used to search for a possible 
spin ordering transition and the Ferrenberg-Swendsen approach was used to cover a range 
of k (and sometimes a) from a single simulation. At a = 3.5, on small lattices (up to 
A^3 = 1600), no peak was observed in the specific heat as a function of the Ising coupling 
and correlation times remained small. Instead, the specific heat was observed to rise to a 
plateau as k was increased. For large k, nearest neighbor spins are strongly correlated and 
the fiuctuation in the number of nodes contributes significantly to the specific heat. This 
is the origin of the plateau. On larger lattices, A''3 = 3200 and A''3 = 6400, a peak appears 
above the plateau value and grows (along with an increased correlation time) suggesting 
the existence of a phase transition. Figure 2 shows the peak just beginning to appear for 
3200 tetrahedra. 

Very long runs (over half a million sweeps) coupled with the Ferrenberg-Swensen 
method allowed us to calculate the kurtosis (1- < > /3 < >^ where E = kj^tt) 



the distribution of the spin correlation for couphngs near that of the specific heat peak [11]. 
Tables 1 and 2 give a sample of the results. The minimum value was near 0.48 for both 
volumes while the value at the peak of the specific heat increased with increasing volume. 
The most naive extrapolation to infinite volume still leaves the kurtosis significantly below 
its gaussian value (2/3) indicating a first order transition. The peak of the specific heat 
itself roughly doubles as the volume is doubled, another hint that the transition is first 
order. 

What happens as a is increased? For significantly greater values of a we were unable 
to find a peak rising above the plateau. It is possible that the first order transition ends 
and is replaced with a smooth crossover behavior, but it appears more likely that it is 
simply necessary to go to larger lattices in order to continue to see the transition. 

In conclusion, the ferromagnetic sector of the phase diagram of three dimensional 
quantum gravity coupled to Ising spins appears to have three components (see figure 1). 
Phase 1 has ordered Ising spins and negative curvature. Phase 2 has ordered Ising spins 
and positive curvature. Phase 3 has disordered Ising spins and positive curvature. Near 
the gravitational transition, the transition from phase 2 to phase 3 appears to be first 
order. For larger a the nature of the transition is uncertain. Further calculations on larger 
lattices are needed to resolve this question and to verify the order of the transition. 

This work was supported, in part, by NSF grant PHY 87-01775. Some calculations 
were performed on the Florida State University Cray Y-MP and at the National Center for 
Supercomputer Applications in Illinois. Also, we acknowledge National Science Foundation 
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FIGURE CAPTIONS 



[1] The phase diagram suggested by quenched mean field theory. 

[2] The specific heat, C, versus the Ising couphng for 3200 tetrahedra. 
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Table 1. Ferrenberg-Swendsen estimates of the specific heat, C, and kurtosis, 
function of Ising coupling for 3200 tetrahedra. The actual coupling for the run was 0.120. 



k 


C 


K, 


0.0790 


0.029 


0.653 


0.0820 


0.062 


0.568 


0.0844 


0.154 


0.474 


0.0866 


0.248 


0.557 


0.0900 


0.124 


0.613 


0.0920 


0.092 


0.608 


0.0960 


0.079 


0.631 



Table 2. Ferrenberg-Swendsen estimates of the specific heat, C, and kurtosis, K/^ clS cL 
function of Ising coupling for 6400 tetrahedra. The actual coupling for the run was 0.090. 
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